home
***
CD-ROM
|
disk
|
FTP
|
other
***
search
/
PD Collection CD 1
/
PD Collection CD 1.iso
/
programer2
/
pari2
/
pari
/
c
/
elliptic
< prev
next >
Wrap
Text File
|
1991-11-28
|
54KB
|
1,769 lines
/********************************************************************/
/********************************************************************/
/** **/
/** COURBES ELLIPTIQUES **/
/** **/
/** Copyright Babe Cool **/
/** **/
/********************************************************************/
/********************************************************************/
#include "genpari.h"
static int aux(),aux2(),numroots2(),numroots3();
GEN smallinitell(x)
GEN x;
{
GEN y,b2,b4,b6,b8,d,j,a11,a13,a33,a64,b81,b22,c4,c6;
long i,av,tetpil;
av=avma;y=cgetg(14,17);
if ((typ(x) != 17) || (lg(x) < 6)) err(elliper1);
for(i=1;i<=5;i++) y[i]=x[i];
b2=gadd(a11=gmul(y[1],y[1]),gmul2n(y[2],2));y[6]=(long)b2;
b4=gadd(a13=gmul(y[1],y[3]),gmul2n(y[4],1));y[7]=(long)b4;
b6=gadd(a33=gmul(y[3],y[3]),a64=gmul2n(y[5],2));y[8]=(long)b6;
b81=gadd(gadd(gmul(a11,y[5]),gmul(a64,y[2])),gmul(y[2],a33));
b8=gsub(b81,gmul(y[4],gadd(y[4],a13)));y[9]=(long)b8;
c4=gadd(b22=gmul(b2,b2),gmulsg(-24,b4));y[10]=(long)c4;
c6=gadd(gmul(b2,gsub(gmulsg(36,b4),b22)),gmulsg(-216,b6));y[11]=(long)c6;
b81=gadd(gmul(b22,b8),gmulsg(27,gmul(b6,b6)));
d=gsub(gmul(b4,gadd(gmulsg(9,gmul(b2,b6)),gmulsg(-8,gmul(b4,b4)))),b81);
y[12]=(long)d;j=gdiv(gmul(gmul(c4,c4),c4),d);y[13]=(long)j;
tetpil=avma;return gerepile(av,tetpil,gcopy(y));
}
static GEN initellcom(x,prec,hard)
GEN x;
long prec;
long hard;
{
GEN y,b2,b4,b6,b8,d,j,a11,a13,a33,a64,b81,b22,c4,c6,p1,p2,p,u,w,w2;
GEN aa0,bb0,aa1,bb1,r1,x0,x1,u2,q,pv,bmod0,bmod1,tau,pi2,pitemp;
long ty,i,av,tetpil,e0,e1,e4,ind,alpha,sw;
av=avma;y=cgetg(20,17);
if ((typ(x) != 17) || (lg(x) < 6)) err(elliper1);
e0=BIGINT;
for(i=1;i<=5;i++)
{
q=(GEN)x[i];y[i]=(long)q;
if(typ(q)==7)
{
e0=min(e0,signe(q[4])?precp(q)+valp(q):valp(q));
p=(GEN)q[2];
}
}
if(e0<BIGINT)
{
q=ggrando(p,e0);for(i=1;i<=5;i++) y[i]=ladd(q,x[i]);
}
b2=gadd(a11=gmul(y[1],y[1]),gmul2n(y[2],2));y[6]=(long)b2;
b4=gadd(a13=gmul(y[1],y[3]),gmul2n(y[4],1));y[7]=(long)b4;
b6=gadd(a33=gmul(y[3],y[3]),a64=gmul2n(y[5],2));y[8]=(long)b6;
b81=gadd(gadd(gmul(a11,y[5]),gmul(a64,y[2])),gmul(y[2],a33));
b8=gsub(b81,gmul(y[4],gadd(y[4],a13)));y[9]=(long)b8;
c4=gadd(b22=gmul(b2,b2),gmulsg(-24,b4));y[10]=(long)c4;
c6=gadd(gmul(b2,gsub(gmulsg(36,b4),b22)),gmulsg(-216,b6));y[11]=(long)c6;
b81=gadd(gmul(b22,b8),gmulsg(27,gmul(b6,b6)));
d=gsub(gmul(b4,gadd(gmulsg(9,gmul(b2,b6)),gmulsg(-8,gmul(b4,b4)))),b81);
y[12]=(long)d;j=gdiv(gmul(gmul(c4,c4),c4),d);y[13]=(long)j;
if(gcmp0(d)) err(initeler2);
ty=typ(d);
if(prec&&(ty<=8)&&(ty!=3)&&(ty!=7))
{
p1=cgetg(6,10);p1[1]=0x1000006;p1[2]=(long)b6;
p1[3]=lmul2n(b4,1);p1[4]=(long)b2;p1[5]=lstoi(4);
if(hard) p1=rootslong(p1,prec);else p1=roots(p1,prec);
if(gsigne(d)>0)
{
p1=greal(p1);
ind=1;e4=p1[1];for(i=2;i<=3;i++)
if(gcmp(p1[i],e4)>0) {ind=i;e4=p1[i];}
p1[ind]=p1[1];p1[1]=e4;
if(gcmp(p1[2],p1[3])<0) {e4=p1[3];p1[3]=p1[2];p1[2]=e4;}
}
else p1[1]=lreal(p1[1]);
y[14]=(long)p1;e1=p1[1];p2=gadd(gmulsg(3,e1),gmul2n(b2,-2));
w=gsqrt(gmul2n(gadd(b4,gmul(e1,gadd(b2,gmulsg(6,e1)))),1),prec);
if(gsigne(p2)>0) w=gneg(w);
aa1=gmul2n(gsub(w,p2),-2);sw=signe(w);
bb1=gmul2n(w,-1);r1=gsub(aa1,bb1);x1=gmul2n(r1,-2);
do
{
aa0=aa1;bb0=bb1;x0=x1;
bb1=gsqrt(gmul(aa0,bb0),prec);setsigne(bb1,sw);
aa1=gmul2n(gadd(gadd(aa0,bb0),gmul2n(bb1,1)),-2);
r1=gsub(aa1,bb1);
x1=gmul(x0,gsqr(gmul2n(gaddsg(1,gsqrt(gdiv(gadd(x0,r1),x0),prec)),-1)));
}
while(gexpo(r1)>=gexpo(bb1)-((prec-2)<<5)+5);
u2=ginv(gmul2n(aa1,2));
w=gaddsg(1,ginv(gmul2n(gmul(u2,x1),1)));
if(gsigne(greal(w))>0) q=ginv(gadd(w,gsqrt(gaddgs(gmul(w,w),-1),prec)));
else q=gsub(w,gsqrt(gaddgs(gmul(w,w),-1),prec));
pitemp=mppi(prec);pi2=gmul2n(pitemp,1);if(gexpo(q)>=0) q=ginv(q);
tau=gmul(gdiv(glog(q,prec),pi2),gneg(gi));
y[19]=lmul(gmul(gsqr(pi2),gabs(u2,prec)),gimag(tau));
u=gmul(pi2,gsqrt(gneg(u2),prec));w2=gmul(tau,u);
if(sw<0) y[15]=(long)u;
else
{
y[15]=lmul2n(gabs(w2[1],prec),1);
q=gexp(gmul(gmul(pi2,gi),gdiv(w2,y[15])),prec);
}
y[16]=(long)w2;
p1=gdiv(gmul(pitemp,pitemp),gmulsg(6,y[15]));q=gsqrt(q,prec);
y[17]=lmul(p1,gdiv(thetanullk(q,3,prec),thetanullk(q,1,prec)));
y[18]=ldiv(gsub(gmul(y[17],y[16]),gmul(gi,pitemp)),y[15]);
}
else
if(ty==7)
{
if(valp(j)>=0) err(initeler1);
p=(GEN)d[2];alpha=valp(c4)>>1;
setvalp(c4,0);setvalp(c6,0);e1=ldivgs(gdiv(c6,c4),6);
c4=gdivgs(c4,48);c6=gdivgs(c6,864);
do
{
e0=e1;p2=gmul(e0,e0);
e1=ldiv(gadd(gmul2n(gmul(e0,p2),1),c6),gsub(gmulsg(3,p2),c4));
}
while(!gegal(e0,e1));
setvalp(e1,valp(e1)+alpha);e1=lsub(e1,gdivgs(b2,12));
w=gsqrt(gmul2n(gadd(b4,gmul(e1,gadd(b2,gmulsg(6,e1)))),1));
p1=gaddgs(gdiv(gmulsg(3,e0),w),1);
if(valp(p1)<=0) w=gneg(w);
y[18]=(long)w;
aa1=gmul2n(gsub(w,gadd(gmulsg(3,e1),gmul2n(b2,-2))),-2);
bb1=gmul2n(w,-1);r1=gsub(aa1,bb1);x1=gmul2n(r1,-2);
pv=gegal(p,deux) ? stoi(4) : p;bmod0=modii(bb1[4],pv);
do
{
aa0=aa1;bb0=bb1;x0=x1;
bb1=gsqrt(gmul(aa0,bb0));bmod1=modii(bb1[4],pv);
if(!gegal(bmod1,bmod0)) bb1=gneg(bb1);
aa1=gmul2n(gadd(gadd(aa0,bb0),gmul2n(bb1,1)),-2);
r1=gsub(aa1,bb1);
p1=gsqrt(gdiv(gadd(x0,r1),x0));
if(!gegal(modii(p1[4],pv),un)) p1=gneg(p1);
x1=gmul(x0,gsqr(gmul2n(gaddsg(1,p1),-1)));
}
while(!gcmp0(r1));
u2=ginv(gmul2n(aa1,2));
w=gaddsg(1,ginv(gmul2n(gmul(u2,x1),1)));
q=ginv(gadd(w,gsqrt(gaddgs(gmul(w,w),-1))));
if(valp(q)<0) q=ginv(q);
y[15]=(long)u2;
y[16]=((kronecker(u2[4],p)>0)&&(!(valp(u2)&1)))?lsqrt(u2):zero;
y[17]=(long)q;p1=cgetg(2,17);
y[14]=(long)p1;p1[1]=e1;y[19]=zero;
}
else {y[14]=y[15]=y[16]=y[17]=y[18]=y[19]=zero;}
tetpil=avma;return gerepile(av,tetpil,gcopy(y));
}
GEN initell2(x,prec)
GEN x;
long prec;
{
return initellcom(x, prec, 0);
}
GEN initell(x,prec)
GEN x;
long prec;
{
return initellcom(x, prec, 1);
}
GEN coordch(e,ch)
GEN e,ch;
{
GEN y,p1,p2,v,v2,v3,v4,v6,r,s,t,u;
long av,tetpil,i,lx=lg(e);
if ((typ(e) != 17) || (lx < 14)) err(elliper1);
u=(GEN)ch[1];r=(GEN)ch[2];s=(GEN)ch[3];t=(GEN)ch[4];
av=avma;y=cgetg(lx,17);v=ginv(u);
y[1]=lmul(v,gadd(e[1],gmul2n(s,1)));v2=gmul(v,v);
y[2]=lmul(v2,gsub(gadd(e[2],gmulsg(3,r)),gmul(s,gadd(e[1],s))));
v3=gmul(v,v2);
y[3]=lmul(v3,p1=gadd(gadd(gmul(r,e[1]),gmul2n(t,1)),e[3]));
v4=gmul(v,v3);
p1=gsub(e[4],gadd(gmul(t,e[1]),gmul(s,p1)));
y[4]=lmul(v4,gadd(p1,gmul(r,gadd(gmul2n(e[2],1),gmulsg(3,r)))));
p1=gadd(e[5],gmul(r,gadd(e[4],gmul(r,gadd(e[2],r)))));
p2=gmul(t,gadd(gadd(gmul(r,e[1]),t),e[3]));v6=gmul(v2,v4);
y[5]=lmul(v6,gsub(p1,p2));
y[6]=lmul(v2,gadd(e[6],gmulsg(12,r)));
y[7]=lmul(v4,gadd(e[7],gmul(r,gadd(e[6],gmulsg(6,r)))));
y[8]=lmul(v6,gadd(e[8],gmul(r,gadd(gmul2n(e[7],1),gmul(r,gadd(e[6],gmul2n(r,2)
))))));
p1=gadd(gmulsg(3,e[7]),gmul(r,gadd(e[6],gmulsg(3,r))));
y[9]=lmul(gmul(v6,v2),gadd(e[9],gmul(r,gadd(gmulsg(3,e[8]),gmul(r,p1)))));
y[10]=lmul(v4,e[10]);y[11]=lmul(v6,e[11]);
y[12]=lmul(gmul(v6,v6),e[12]);y[13]=lcopy(e[13]);
if(lx==14) {tetpil=avma;return gerepile(av,tetpil,gcopy(y));}
p1=(GEN)e[14];
if (gcmp0(p1))
{
y[14] = y[15] = y[16] = y[17] = y[18] = y[19] = zero;
}
else
{
p2=cgetg(4,18);
for(i=1;i<=3;i++) p2[i]=lmul(v2,gsub(p1[i],r));
y[14]=(long)p2;
y[15]=lmul(u,e[15]);if(typ(e[1])==7) y[15]=lmul(u,y[15]);
y[16]=lmul(u,e[16]);
y[17]=ldiv(e[17],u);y[18]=ldiv(e[18],u);
y[19]=lmul(gmul(u,u),e[19]);
}
tetpil=avma;return gerepile(av,tetpil,gcopy(y));
}
GEN pointch(x,ch)
GEN x,ch;
{
GEN y,p1,p2,p3,v,v2,v3,mor,r,s,t,u;
long av,tetpil,tx,lx=lg(x),i;
if ((typ(x) != 17) || (typ(ch) != 17)) err(elliper1);
if (lx < 2) return gcopy(x);
av=avma;u=(GEN)ch[1];r=(GEN)ch[2];s=(GEN)ch[3];t=(GEN)ch[4];
tx=typ(x[1]);v=ginv(u);v2=gmul(v,v);v3=gmul(v,v2);
mor=gneg(r);
if(tx>=17)
{
y=cgetg(lx,tx);
for(i=1;i<lx;i++)
{
p2=(GEN)x[i];
if(lg(p2)<3) y[i]=(long)p2;
else
{
p1=cgetg(3,17);
p1[1]=lmul(v2,p3=gadd(p2[1],mor));
p1[2]=lmul(v3,gsub(p2[2],gadd(gmul(s,p3),t)));
y[i]=(long)p1;
}
}
}
else
{
if(lg(x)<3) y=x;
else
{
y=cgetg(3,17);
y[1]=lmul(v2,p3=gadd(x[1],mor));
y[2]=lmul(v3,gsub(x[2],gadd(gmul(s,p3),t)));
}
}
tetpil=avma;return gerepile(av,tetpil,gcopy(y));
}
int oncurve(e,z)
GEN e,z;
{
GEN p1,p2,x,y;
long av=avma,f;
if ((typ(e) != 17) || (lg(e) < 6) || (typ(z) != 17)) err(elliper1);
if(lg(z)<3) return 1; /* cas du point a l'infini */
x=(GEN)z[1];y=(GEN)z[2];
p1=gmul(y,gadd(y,gadd(e[3],gmul(e[1],x))));
p2=gadd(e[5],gmul(x,gadd(e[4],gmul(x,gadd(e[2],x)))));
p1=gsub(p1,p2);
f=precision(p1);
if(!f) {f=gcmp0(p1);avma=av;return f;}
f=gcmp0(p1)||(gexpo(p1)< -f+5);avma=av;return f;
}
GEN hells(e,x,prec)
GEN e,x;
long prec;
{
GEN w,z,t,mu,unreel;
long av=avma,tetpil,n,lim;
if(!oncurve(e,x)) err(heller1);
affsr(1,unreel=cgetr(prec));
t=gdiv(unreel,x[1]);mu=gmul2n(glog(numer(x[1]),prec),-1);
n=0;lim=(prec-2)<<4;
do
{
w=gmul(t,gaddsg(4,gmul(t,gadd(e[6],gmul(t,gadd(gmul2n(e[7],1),gmul(t,e[8])))))
));
z=gsub(gun,gmul(gmul(t,t),gadd(e[7],gmul(t,gadd(gmul2n(e[8],1),gmul(t,e[9]))))
));
mu=gadd(mu,gmul2n(glog(z,prec),-((n<<1)+3)));t=gdiv(w,z);n++;
}
while(n<=(lim+5));
tetpil=avma;return gerepile(av,tetpil,gcopy(mu));
}
GEN hell2(e,x,prec)
GEN e,x;
long prec;
{
GEN ep,e3,ro,p1,p2,mu,d,xp;
long av=avma,tetpil,lx,lc,i,j,tx;
if((typ(e)!=17)||(lg(e)<20)) err(elliper1);
if(!oncurve(e,x)) err(heller1);
d=(GEN)e[12];ro=(GEN)e[14];e3=(gsigne(d)<0)?(GEN)ro[1]:(GEN)ro[3];
p1=cgetg(5,17);p1[1]=un;p1[2]=laddgs(gfloor(e3),-1);
p1[3]=p1[4]=zero;ep=coordch(e,p1);xp=pointch(x,p1);
if(typ(x[1])<17)
{
if(lg(x)<3) return gzero;
tetpil=avma;return gerepile(av,tetpil,hells(ep,xp,prec));
}
else
{
lx=lg(x);tetpil=avma;mu=cgetg(lx,tx=typ(x));
if(tx<19) for(i=1;i<lx;i++) mu[i]=(long)hells(ep,xp[i],prec);
else
{
lc=lg(x[1]);
for(i=1;i<lx;i++)
{
p1=cgetg(lc,18);mu[i]=(long)p1;p2=(GEN)xp[i];
for(j=1;j<lc;j++) p1[j]=(long)hells(ep,p2[j],prec);
}
}
return gerepile(av,tetpil,mu);
}
}
GEN addell(e,z1,z2)
GEN e,z1,z2;
{
GEN y,p1,p2,x1,x2,x3,y1,y2,y3,al;
long av=avma,tetpil;
if((typ(e) != 17) || (typ(z1) != 17) || (typ(z2) != 17)) err(elliper1);
if(lg(z1)<3) return gcopy(z2);
if(lg(z2)<3) return gcopy(z1);
x1=(GEN)z1[1];x2=(GEN)z2[1];y1=(GEN)z1[2];y2=(GEN)z2[2];
if(gegal(x1,x2))
if(!gegal(y1,y2)) {y=cgetg(2,17);y[1]=zero;return y;}
else
{
p1=gadd(gsub(e[4],gmul(e[1],y1)),gmul(x1,gadd(gmul2n(e[2],1),gmulsg(3,x1))));
p2=gadd(e[3],gadd(gmul(e[1],x1),gmul2n(y1,1)));
if(gcmp0(p2)) {avma=av;y=cgetg(2,17);y[1]=zero;return y;}
}
else {p1=gsub(y2,y1);p2=gsub(x2,x1);}
al=gdiv(p1,p2);
x3=gsub(gmul(al,gadd(al,e[1])),gadd(gadd(x1,x2),e[2]));
y3=gneg(gadd(gadd(gadd(e[3],y1),gmul(x3,e[1])),gmul(al,gsub(x3,x1))));
tetpil=avma;y=cgetg(3,17);y[1]=lcopy(x3);y[2]=lcopy(y3);
return gerepile(av,tetpil,y);
}
GEN subell(e,z1,z2)
GEN e,z1,z2;
{
GEN zp;
long av=avma,tetpil;
if((typ(e) != 17) || (typ(z2) != 17)) err(elliper1);
if(lg(z2)<3) return gcopy(z1);
zp=cgetg(3,17);zp[1]=z2[1];
zp[2]=lneg(gadd(gadd(z2[2],e[3]),gmul(e[1],z2[1])));
tetpil=avma;return gerepile(av,tetpil,addell(e,z1,zp));
}
GEN ordell(e,x,prec)
GEN e,x;
long prec;
{
GEN p1,p2,p3,d,y;
long av=avma,tetpil,td,i,lx,tx=typ(x);
if((typ(e) != 17) || (lg(e) < 6)) err(elliper1);
if(tx<17)
{
p1=gadd(e[5],gmul(x,gadd(e[4],gmul(x,gadd(e[2],x)))));
p2=gadd(e[3],gmul(e[1],x));d=gadd(gmul(p2,p2),gmul2n(p1,2));
if((td=typ(d))==1)
{
if(gsigne(d)<0) {avma=av;return cgetg(1,17);}
p3=racine(d);
if(!gegal(d,gmul(p3,p3))) {avma=av;return cgetg(1,17);}
p1=gsub(p3,p2);tetpil=avma;
if(signe(d)) {y=cgetg(3,17);y[1]=lmul2n(p1,-1);y[2]=lsub(y[1],p3);}
else {y=cgetg(2,17);y[1]=lmul2n(p1,-1);}
return gerepile(av,tetpil,y);
}
else
if((td==3)&&gegal(d[1],gdeux))
{
if(gcmp0(p2))
{
avma=av;y=cgetg(2,17);p3=cgetg(3,3);y[1]=(long)p3;
p3[1]=deux;p3[2]=gcmp0(p1)?zero:un;
}
else
if(gcmp0(p1))
{
avma=av;y=cgetg(3,17);p3=cgetg(3,3);y[1]=(long)p3;
p3[1]=deux;p3[2]=zero;p3=cgetg(3,3);y[2]=(long)p3;
p3[1]=deux;p3[2]=un;
}
else {avma=av;y=cgetg(1,17);}
return y;
}
else
{
if((td==3)&&(kronecker(d[2],d[1])== -1))
{avma=av;y=cgetg(1,17);return y;}
p3=gsqrt(d,prec);p1=gsub(p3,p2);tetpil=avma;
if(!gcmp0(d))
{
y=cgetg(3,17);y[1]=lmul2n(p1,-1);y[2]=lsub(y[1],p3);
}
else {y=cgetg(2,17);y[1]=lmul2n(p1,-1);}
return gerepile(av,tetpil,y);
}
}
else
{
lx=lg(x);y=cgetg(lx,tx);
for(i=1;i<lx;i++) y[i]=(long)ordell(e,x[i],prec);
return y;
}
}
GEN powell(e,n,z)
GEN e,n,z;
{
GEN y,zp;
long s=signe(n),av=avma,i,j,tetpil;
unsigned long m;
if(typ(e) != 17) err(elliper1);
if(!s) {y=cgetg(2,17);y[1]=zero;return y;}
if(lg(z)<3) return gcopy(z);
if(s<0)
{
n=gneg(n);zp=cgetg(3,17);zp[1]=z[1];
zp[2]=lneg(gadd(gadd(z[2],e[3]),gmul(e[1],z[1])));
}
else zp=z;
if(gcmp1(n)) {tetpil=avma;return gerepile(av,tetpil,gcopy(zp));}
y=cgetg(2,17);y[1]=zero;
for (i=lgef(n)-1;i>2;i--)
for (m=n[i],j=0;j<32;j++,m>>=1)
{
if (m&1) y=addell(e,y,zp);
zp=addell(e,zp,zp);
}
for (m=n[2];m>1;m>>=1)
{
if (m&1) y=addell(e,y,zp);
zp=addell(e,zp,zp);
}
tetpil=avma;y=addell(e,y,zp);
return gerepile(av,tetpil,y);
}
GEN matell(e,x,prec)
GEN e,x;
long prec;
{
GEN y,p1,p2,pdiag;
long av=avma,tetpil,lx=lg(x),i,j;
if(typ(x) != 17) err(elliper1);
lx=lg(x);y=cgetg(lx,19);pdiag=cgetg(lx,18);
for(i=1;i<lx;i++) {pdiag[i]=(long)ghell(e,x[i],prec);y[i]=lgetg(lx,18);}
for(i=1;i<lx;i++)
{
p1=(GEN)y[i];p1[i]=lmul2n(pdiag[i],1);
for(j=i+1;j<lx;j++)
{
p2=gsub(ghell(e,addell(e,x[i],x[j]),prec),gadd(pdiag[i],pdiag[j]));
p1[j]=(long)p2;((GEN)y[j])[i]=(long)p2;
}
}
tetpil=avma;return gerepile(av,tetpil,gcopy(y));
}
GEN zell(e,z,prec)
GEN e,z;
long prec;
{
long av=avma,tetpil,ty,sw,fl;
GEN t,u,w,p1,p2,b2,b4,r0,r1,aa1,aa0,bb1,bb0,x0,x1,c0,e1,delta;
GEN bmod0,bmod1,p,u2;
if((typ(e)!=17)||(lg(e)<20)) err(elliper1);
if(!oncurve(e,z)) err(heller1);
ty=typ(e[12]);b2=(GEN)e[6];b4=(GEN)e[7];
if(lg(z)<3) return (ty==7)?gun:gzero;
e1=(GEN)((GEN)e[14])[1];
w=gsqrt(gmul2n(gadd(b4,gmul(e1,gadd(b2,gmulsg(6,e1)))),1),prec);
if((ty<=8)&&(ty!=3)&&(ty!=7))
{
e1=(GEN)((GEN)e[14])[1];p2=gadd(gmulsg(3,e1),r0=gmul2n(b2,-2));
if(gsigne(p2)>0) w=gneg(w);
aa1=gmul2n(gsub(w,p2),-2);sw=signe(w);r0=gadd(e1,r0);
bb1=gmul2n(w,-1);r1=gsub(aa1,bb1);
c0=gadd(z[1],gmul2n(r0,-1));
if(gsigne(c0))
{
delta=gdiv(gmul(aa1,r1),gsqr(c0));
x0=gmul2n(gmul(c0,gaddsg(1,gsqrt(gsubsg(1,gmul2n(delta,2)),prec))),-1);
}
else x0=gsqrt(gneg(gmul(aa1,r1)),prec);
x1=gmul(x0,gsqr(gmul2n(gaddsg(1,gsqrt(gdiv(gadd(x0,r1),x0),prec)),-1)));
fl=0;
do
{
aa0=aa1;bb0=bb1;x0=x1;r0=r1;
bb1=gsqrt(gmul(aa0,bb0),prec);setsigne(bb1,sw);
aa1=gmul2n(gadd(gadd(aa0,bb0),gmul2n(bb1,1)),-2);
r1=gsub(aa1,bb1);
x1=gmul(x0,gsqr(gmul2n(gaddsg(1,gsqrt(gdiv(gadd(x0,r1),x0),prec)),-1)));
if(gexpo(gsub(x1,x0))<gexpo(x1)-((prec-2)<<5)+5) fl++;
}
while(fl<2);
t=gaddsg(1,u=gdiv(x1,aa1));
if(gcmp0(t)||(gexpo(t)<-((prec-2)<<5)+5)) t=gneg(gun);
else t=gdiv(u,gsqr(gaddsg(1,gsqrt(t,prec))));
u=gsqrt(ginv(gmul2n(aa1,2)),prec);t=glog(t,prec);t=gmul(t,u);
/* les deux lignes suivantes ont ete rajoutees au pif. A verifier */
if(gsigne(c0)*gsigne(gadd(e[3],gadd(gmul(e[1],z[1]),gmul2n(z[2],1))))<0)
t=gneg(t);
p1=gsub(p2=gdiv(gimag(t),((GEN)e[16])[2]),gmul2n(gun,-2));
if(gcmp(gabs(p1,prec),ghalf)>=0)
{
t=gsub(t,gmul(e[16],gfloor(gadd(p2,lisexpr("0.1")))));
}
if(signe(greal(t))<0) t=gadd(t,e[15]);
tetpil=avma;return gerepile(av,tetpil,gcopy(t));
}
else
if(ty==7)
{
w=(GEN)e[18];r0=gadd(e1,gmul2n(b2,-2));p=(GEN)((GEN)e[12])[2];
aa1=gmul2n(gsub(w,gadd(gmulsg(3,e1),gmul2n(b2,-2))),-2);
bb1=gmul2n(w,-1);r1=gsub(aa1,bb1);
c0=gadd(z[1],gmul2n(r0,-1));delta=gdiv(gmul(aa1,r1),gsqr(c0));
x0=gmul2n(gmul(c0,gaddsg(1,gsqrt(gsubsg(1,gmul2n(delta,2)),prec))),-1);
bmod0=modii(bb1[4],p);
x1=gmul(x0,gsqr(gmul2n(gaddsg(1,gsqrt(gdiv(gadd(x0,r1),x0),prec)),-1)));
do
{
aa0=aa1;bb0=bb1;x0=x1;r0=r1;
bb1=gsqrt(gmul(aa0,bb0));bmod1=modii(bb1[4],p);
if(!gegal(bmod1,bmod0)) bb1=gneg(bb1);
aa1=gmul2n(gadd(gadd(aa0,bb0),gmul2n(bb1,1)),-2);
r1=gsub(aa1,bb1);
p1=gsqrt(gdiv(gadd(x0,r1),x0));
if(!gegal(modii(p1[4],p),un)) p1=gneg(p1);
x1=gmul(x0,gsqr(gmul2n(gaddsg(1,p1),-1)));
}
while(!gcmp0(r1));
u2=ginv(gmul2n(aa1,2));
if(!gcmp0(e[16]))
{
t=gsqrt(gaddsg(1,gdiv(x1,aa1)),prec);
t=gdiv(gaddsg(-1,t),gaddsg(1,t));
}
else t=gaddsg(2,ginv(gmul(u2,x1)));
tetpil=avma;return gerepile(av,tetpil,gcopy(t));
}
else err(zeller1);
}
GEN pointell(e,z,prec)
GEN e,z;
long prec;
{
long av=avma,tetpil;
GEN p1,pii2,pii4,pii6,a,tau,q,u,y,yp,u1,u2,qn,qnu2,qnu1,qnu,qnu4,qnu3,p2,v;
if((typ(e)!=17)||(lg(e)<20)) err(elliper1);
p1=mppi(prec);setexpo(p1,2);
pii2=cgetg(3,6);pii2[1]=zero;pii2[2]=(long)p1;
z=gdiv(z,e[15]);tau=gdiv(e[16],e[15]);
a=ground(gdiv(gimag(z),gimag(tau)));z=gsub(z,gmul(a,tau));
a=ground(greal(z));z=gsub(z,a);
if(gcmp0(z)||gexpo(gnorm(z))< -((prec-2)<<6)+10) {avma=av;v=cgetg(2,17);v[1]=zero;return v;}
q=gexp(gmul(pii2,tau),prec);u=gexp(gmul(pii2,z),prec);
y=gadd(gdivgs(gun,12),gdiv(u,u2=gsqr(u1=gsubsg(1,u))));
yp=gdiv(gaddsg(1,u),gmul(u1,u2));
qn=q;pii2=gdiv(pii2,e[15]);pii4=gmul(pii2,pii2);pii6=gmul(pii4,pii2);
do
{
p1=gsub(gmul(u,gadd(gdivsg(1,qnu2=gsqr(qnu1=gsubsg(1,qnu=gmul(qn,u)))),gdivsg(1,qnu4=gsqr(qnu3=gsub(qn,u))))),gmul2n(gdivsg(1,gsqr(gsubsg(1,qn))),1));
y=gadd(y,gmul(qn,p1));
p2=gadd(gdiv(gaddsg(1,qnu),gmul(qnu1,qnu2)),gdiv(gadd(qn,u),gmul(qnu3,qnu4)));
yp=gadd(yp,gmul(qn,p2));qn=gmul(q,qn);
}
while(gexpo(qn)>-((prec-2)<<5)-5);
yp=gmul(u,gmul(pii6,yp));y=gsub(gmul(pii4,y),gdivgs(e[6],12));
yp=gsub(yp,gadd(e[3],gmul(e[1],y)));
tetpil=avma;v=cgetg(3,17);v[1]=lcopy(y);v[2]=lmul2n(yp,-1);
return gerepile(av,tetpil,v);
}
GEN apell2(e,p)
GEN e,p;
{
/* Calcul de a_p par les sommes de Jacobi */
GEN y,e71;
long av=avma,av2,i,l,s,e72,e6,e8;
if((typ(e)!=17)||(lg(e)<14)) err(elliper1);
if(lgef(p)>3) err(apeller1);
s=0;l=p[2];e71=gmul2n(e[7],1);av2=avma;
if((unsigned long)l>=0x40000000) err(apeller1);
if((l<757)&&(l>2))
{
e6=itos(modis(e[6],l));e72=(itos(modis(e[7],l)))<<1;
e8=itos(modis(e[8],l));
for(i=0;i<l;i++) s=s+kross(e8+i*(e72+i*(e6+(i<<2))),l);
}
else
{
for(i=0;i<l;i++)
{
y=gadd(e[8],gmulsg(i,gadd(e71,gmulsg(i,gaddsg(i<<2,e[6])))));
s=s+krogs(y,l);avma=av2;
}
}
avma=av;return stoi(-s);
}
GEN addsell(e,z1,z2)
GEN e,z1,z2;
{
GEN y,p1,p2,x1,x2,x3,y1,y2,y3,al;
long av=avma,tetpil;
if(lg(z1)<3) return gcopy(z2);
if(lg(z2)<3) return gcopy(z1);
x1=(GEN)z1[1];x2=(GEN)z2[1];y1=(GEN)z1[2];y2=(GEN)z2[2];
if(gegal(x1,x2))
{
if(!gegal(y1,y2)) {y=cgetg(2,17);y[1]=zero;return y;}
else
{
p1=gadd(e,gmul(x1,gmulsg(3,x1)));
p2=gmul2n(y1,1);
if(gcmp0(p2)) {avma=av;y=cgetg(2,17);y[1]=zero;return y;}
}
}
else {p1=gsub(y2,y1);p2=gsub(x2,x1);}
al=gdiv(p1,p2);
x3=gsub(gmul(al,al),gadd(x1,x2));
y3=gneg(gadd(y1,gmul(al,gsub(x3,x1))));
tetpil=avma;y=cgetg(3,17);y[1]=lcopy(x3);y[2]=lcopy(y3);
return gerepile(av,tetpil,y);
}
GEN doubsell(e,z1)
GEN e,z1;
{
GEN x1,x3,y3,y,y1,p1,p2,al;
long av=avma,tetpil;
if(lg(z1)<3) return gcopy(z1);
x1=(GEN)z1[1];y1=(GEN)z1[2];
p1=gadd(e,gmul(x1,gmulsg(3,x1)));p2=gmul2n(y1,1);
if(gcmp0(p2)) {avma=av;y=cgetg(2,17);y[1]=zero;return y;}
al=gdiv(p1,p2);
x3=gsub(gmul(al,al),gadd(x1,x1));
y3=gneg(gadd(y1,gmul(al,gsub(x3,x1))));
tetpil=avma;y=cgetg(3,17);y[1]=lcopy(x3);y[2]=lcopy(y3);
return gerepile(av,tetpil,y);
}
GEN subsell(e,z1,z2)
GEN e,z1,z2;
{
GEN zp;
long av=avma,tetpil;
if(lg(z2)<3) return gcopy(z1);
zp=cgetg(3,17);zp[1]=z2[1];
zp[2]=lneg(z2[2]);
tetpil=avma;return gerepile(av,tetpil,addsell(e,z1,zp));
}
GEN powsell(e,n,z)
GEN e,n,z;
{
GEN y,zp;
long s=signe(n),av=avma,i,j,tetpil;
unsigned long m;
if(!s) {y=cgetg(2,17);y[1]=zero;return y;}
if(lg(z)<3) return gcopy(z);
if(s<0)
{
n=gneg(n);zp=cgetg(3,17);zp[1]=z[1];
zp[2]=lneg(z[2]);
}
else zp=z;
if(gcmp1(n)) {tetpil=avma;return gerepile(av,tetpil,gcopy(zp));}
y=cgetg(2,17);y[1]=zero;
for (i=lgef(n)-1;i>2;i--)
{
for (m=n[i],j=0;j<32;j++,m>>=1)
{
if (m&1) y=addsell(e,y,zp);
zp=doubsell(e,zp);
}
}
for (m=n[2];m>1;m>>=1)
{
if (m&1) y=addsell(e,y,zp);
zp=doubsell(e,zp);
}
tetpil=avma;y=addsell(e,y,zp);
return gerepile(av,tetpil,y);
}
GEN powssell(e,n,z)
GEN e,z;
long n;
{
GEN y,zp;
long av=avma,tetpil;
unsigned long m;
if(!n) {y=cgetg(2,17);y[1]=zero;return y;}
if(lg(z)<3) return gcopy(z);
if(n<0)
{
n= -n;zp=cgetg(3,17);zp[1]=z[1];
zp[2]=lneg(z[2]);
}
else zp=z;
if(n==1) {tetpil=avma;return gerepile(av,tetpil,gcopy(zp));}
y=cgetg(2,17);y[1]=zero;
for (m=n;m>1;m>>=1)
{
if (m&1) y=addsell(e,y,zp);
zp=doubsell(e,zp);
}
tetpil=avma;y=addsell(e,y,zp);
return gerepile(av,tetpil,y);
}
#define HASHSP 255
GEN apell1(e,p)
GEN e,p;
{
long av,av3,tetpil,k,k2,i,j,j1,lim,limite,succes,com,j2,s;
GEN tabla, tablb, count, index, hash;
GEN p1,p2,p3,q,h,hp,f,fh,fg,ftest;
GEN unmodp,pordmin,u,p1p,p2p,acon,bcon,xp,yp,c4,c6,cp4;
long flc,flcc,sucfin,fll,flf,x,pfinal;
if((typ(e)!=17)||(lg(e)<14)) err(elliper1);
if(gcmpgs(p,20)<0) return apell2(e,p);
tabla = newbloc(1000000);
tablb = newbloc(1000000);
count = newbloc(256);
index = newbloc(257);
hash = newbloc(1000000);
av=avma;limite=(av+bot)>>1;
unmodp=gmodulcp(gun,p);c4=gdivgs(gmul(unmodp,e[10]),-48);
c6=gdivgs(gmul(unmodp,e[11]),-864);
pordmin=gceil(gmul2n(gsqrt(p,4),2));p2p=gmul2n(p1p=gaddsg(1,p),1);
x=0;flcc=0;flc=kronecker(c6[2],p);u=c6;acon=gzero;bcon=gun;
sucfin=1;h=p1p;
while(sucfin)
{
while((flc==flcc)||(!flc))
{
x++;u=gadd(c6,gmulsg(x,gaddgs(c4,x*x)));
flc=kronecker(u[2],p);
}
flcc=flc;
s=itos(gceil(gsqrt(gdiv(pordmin,bcon),3)))>>1;
cp4=gmul(c4,yp=gsqr(u));
xp=gmulsg(x,u);f=cgetg(3,17);f[1]=(long)xp;f[2]=(long)yp;
fh=powsell(cp4,h,f);
if (bcon != gun) f=powsell(cp4,bcon,f); /* sic */
p1=fh;flf=1;
for(i=0;i<=HASHSP;i++) count[i]=0;
for(i=0;(i<=s-1)&&flf;i++)
{
if(lg(p1)==3)
{
p2=(GEN)((GEN)p1[1])[2];tabla[i]=p2[lgef(p2)-1];j=tabla[i]&HASHSP;count[j]++;
/* printf("%d ",j);fflush(stdout); */
p2=(GEN)((GEN)p1[2])[2];tablb[i]=p2[lgef(p2)-1];
p1=addsell(cp4,p1,f);
}
else flf=0;
}
/* printf("\nsj:\n"); */
if(flf)
{
fg=powssell(cp4,s,f);ftest=fg;
index[0]=0;for(i=0;i<=HASHSP-1;i++) index[i+1]=index[i]+count[i];
for(i=0;i<=s-1;i++) hash[index[tabla[i]&HASHSP]++]=i;
index[0]=0;for(i=0;i<=HASHSP;i++) index[i+1]=index[i]+count[i];
succes=0;com=1;av3=avma;pfinal=p[lgef(p)-1];
while(!succes)
{
p1=(GEN)((GEN)ftest[1])[2];k=p1[lgef(p1)-1];j=k&HASHSP;
/* printf("%d ",j);fflush(stdout); */
fll=(lg(ftest)==3);
for(j1=index[j];(j1<index[j+1])&&(!succes);j1++)
{
j2=hash[j1];
if((tabla[j2]==k)&&fll)
{
p2=(GEN)((GEN)ftest[2])[2];k2=p2[lgef(p2)-1];
if((tablb[j2]==k2)||(tablb[j2]==pfinal-k2))
{
p1=addsell(cp4,powssell(cp4,j2,f),fh);
succes=gegal(p1[1],ftest[1]);
}
}
}
if(!succes)
{
com++;
if(avma>=limite) ftest=addsell(cp4,ftest,fg);
else
{
tetpil=avma;ftest=gerepile(av3,tetpil,addsell(cp4,ftest,fg));
}
if (lg(ftest)<3) err(apeller2);
}
}
h=addii(h,mulsi(j2,bcon));p2=mulsi(s,mulsi(com,bcon));
h=(!cmpii(((GEN)p1[2])[2],((GEN)ftest[2])[2]))?subii(h,p2):addii(h,p2);
}
else h=addii(h,mulsi(i-1,bcon));
p2=factor(h);
p1=(GEN)p2[1];
p2=(GEN)p2[2];
for(i=1;i<lg(p1);i++)
{
p3=divii(h,p1[i]);fh=powsell(cp4,p3,f);lim=itos(p2[i]);
for(j=1;(j<=lim)&&(lg(fh)<3);j++)
{
h=p3;if(j<lim) {p3=divii(h,p1[i]);fh=powsell(cp4,p3,f);}
}
}
p1=gmodulcp(acon,bcon);p2=gmodulcp(gzero,h);
p1=chinois(p1,p2);acon=(GEN)p1[2];bcon=(GEN)p1[1];
if(gcmp(bcon,pordmin)>=0)
{
q=ground(gdiv(gsub(p1p,acon),bcon));sucfin=0;
hp=addii(mulii(q,bcon),acon);tetpil=avma;
}
else
{
acon=modii(subii(p2p,acon),bcon);
p1=subii(acon,bcon);if(signe(addii(acon,p1))>0) acon=p1;
q=ground(gdiv(gsub(p1p,acon),bcon));
h=addii(mulii(q,bcon),acon);
}
}
killbloc(tabla); killbloc(tablb); killbloc(count);
killbloc(index); killbloc(hash);
return
(flc==1)?gerepile(av,tetpil,gsub(p1p,hp)):gerepile(av,tetpil,gsub(hp,p1p));
}
typedef struct {int isnull; long x; long y;} sellpt;
long addmod(a, b, p)
long a, b, p;
{
long res = a + b;
return (res >= p) ? res - p : res;
}
long submod(a, b, p)
long a, b, p;
{
long res = a - b;
return (res >= 0) ? res : res + p;
}
extern long mulmodll();
#define mulmod(a,b,c) mulmodll(a,b,c)
static long divmod(a, b, p)
long a, b, p;
{
long v1 = 0, v2 = 1, v3, r, oldp = p;
while (b > 1)
{
v3 = v1 - (p / b) * v2; v1 = v2; v2 = v3;
r = p % b; p = b; b = r;
}
if (v2 < 0) v2 += oldp;
return mulmod(a, v2, oldp);
}
void addsell1(e, p, p1, p2, p3)
long e, p;
sellpt *p1, *p2, *p3;
{
long num, den, lambda;
if (p1->isnull) {*p3 = *p2; return;}
if (p2->isnull) {*p3 = *p1; return;}
if (p1->x == p2->x)
if ((p1->y)&&(p1->y == p2->y))
{
num = addmod(e, mulmod(3, mulmod(p1->x, p1->x, p), p), p);
den = addmod(p1->y, p1->y, p);
}
else
{
p3->isnull = 1;
return;
}
else
{
num = submod(p1->y, p2->y, p);
den = submod(p1->x, p2->x, p);
}
lambda = divmod(num, den, p);
p3->x = submod(mulmod(lambda, lambda, p), addmod(p1->x, p2->x, p), p);
p3->y = submod(mulmod(lambda, submod(p2->x, p3->x, p), p), p2->y, p);
p3->isnull = 0;
}
void powssell1(e, p, n, p1, p2)
long e, p, n;
sellpt *p1, *p2;
{
sellpt p4, p3;
p3 = *p1;
if (n < 0) {n = -n; if (p3.y) p3.y = p - p3.y;}
p2->isnull = 1;
for(;;)
{
if (n&1) addsell1(e, p, p2, &p3, p2);
n>>=1;
if (!n) return;
addsell1(e, p, &p3, &p3, &p4);
p3 = p4;
}
}
typedef struct {long x; long y; long i;} multiple;
int compare_multiples(a, b)
multiple *a, *b;
{
return a->x - b->x;
}
GEN apell(e,pl)
GEN e,pl;
{
long av = avma,i,j,com,s;
GEN p1,p2,unmodp;
sellpt f,fh,fg,ftest,f2;
long pordmin,u,p1p,p2p,acon,bcon,xp,yp,c4,c6,cp4;
long flb,flc,flcc,x, q, h, p3, p, l, r, m;
multiple *table;
if((typ(e)!=17)||(lg(e)<14)) err(elliper1);
if (gcmpgs(pl, 0x3fffffff) > 0) return apell1(e, pl);
p = itos(pl);
if (p < 100) return apell2(e, pl);
unmodp = gmodulcp(gun, pl);
c4 = itos(gdivgs(gmul(unmodp, e[10]), -48)[2]);
c6 = itos(gdivgs(gmul(unmodp, e[11]), -864)[2]);
pordmin = 1 + 4 * sqrt((float)p);
p1p = p + 1; p2p = p1p << 1;
x = 0; flcc = 0; flc = kross(c6, p); u = c6; acon = 0; bcon = 1;
h = p1p;
for(;;)
{
while((flc == flcc) || (!flc))
{
x++;
u = addmod(c6, mulmod(x, c4 + mulmod(x, x, p), p), p);
flc = kross(u,p);
}
flcc = flc;
s = sqrt(((float)pordmin)/bcon) / 2;
if (!s) s = 1;
if(bcon==1) table=(multiple *)malloc((s+1)*sizeof(multiple));
yp = mulmod(u, u, p);
cp4 = mulmod(c4, yp, p);
xp = mulmod(x, u, p);
f.isnull = 0; f.x = xp; f.y = yp;
powssell1(cp4, p, h, &f, &fh);
if (bcon > 1) powssell1(cp4, p, bcon, &f, &f2);else f2=f;
for(i=0; i < s; i++)
if (fh.isnull)
{
h += bcon * i;
goto trouve;
}
else
{
table[i].x = fh.x;
table[i].y = fh.y;
table[i].i = i;
addsell1(cp4, p, &fh, &f2, &fh);
}
qsort(table, s, sizeof(multiple), compare_multiples);
powssell1(cp4, p, s, &f2, &fg); ftest = fg;
for(com = 1;; com++)
{
if (ftest.isnull) err(apeller3);
l = 0; r = s;
while (l < r)
{
m = (l + r) >> 1;
if (table[m].x < ftest.x) l = m + 1; else r = m;
}
if ((r < s) && (table[r].x == ftest.x)) break;
addsell1(cp4, p, &ftest, &fg, &ftest);
}
h += table[r].i * bcon;
if (table[r].y == ftest.y) h -= s * com * bcon; else h += s * com *
bcon;
trouve:
p2=factor(stoi(h));
p1=(GEN)p2[1];
p2=(GEN)p2[2];
for(i=1; i < lg(p1); i++)
for(j = ((GEN)p2[i])[2]; j; j--)
{
p3 = h / ((GEN)p1[i])[2];
powssell1(cp4, p, p3, &f, &fh);
if (!fh.isnull) break;
h = p3;
}
flb=0;
if (bcon > 1)
{
p1 = gmodulcp(stoi(acon), stoi(bcon)); p2=gmodulcp(gzero, stoi(h));
p1=chinois(p1,p2);acon=itos((GEN)p1[2]);bcon=((GEN)p1[1])[2];
if((bcon<0)||(lgef(p1[1])>3)) flb=1;
}
else
bcon = h;
if(flb||(bcon >= pordmin))
{
if(flb) h=acon;
else
{
q = ((unsigned long)(p2p + bcon - (acon << 1))) / (bcon << 1);
h = acon + q * bcon;
}
break;
}
else
{
acon = (p2p - acon) % bcon;
if ((acon << 1) > bcon) acon -= bcon;
q = ((unsigned long)(p2p + bcon - (acon << 1))) / (bcon << 1);
h = acon + q * bcon;
}
}
avma = av;free(table);
return stoi((flc == 1) ? p1p - h : h - p1p);
}
GEN anell(e, n)
GEN e;
long n;
{
long p, pk, i, m, av, tetpil, pl[3];
GEN p1, p2, ap, an;
if((typ(e)!=17)||(lg(e)<14)) err(elliper1);
if (n <= 0) return cgetg(1, 17);
an = cgetg(n+1, 17);
an[1] = un;
for(i=2; i <= n; i++) an[i] = 0;
pl[0] = 0x1000003; pl[1] = 0x1010003;
for (p = 2; p <= n; p++) if (!an[p])
{
pl[2] = p;
if (divise(e[12], pl)) /* mauvaise reduction */
switch (((((p & 3) + 1) & 2) -1) * krogs(e[11], p)) /* renvoie (-c6 /
p) */
{
case -1: /* non deployee */
for(m = p; m <= n; m += p) if (an[m / p]) an[m] = lneg(an[m /
p]);
continue;
case 0: /* additive */
for(m = p; m <= n; m += p) an[m] = zero;
continue;
case 1: /* deployee */
for(m = p; m <= n; m += p) if (an[m / p]) an[m] = lcopy(an[m /
p]);
}
else /* bonne reduction */
for(pk = p; pk <= n; pk *= p)
{
if (pk == p)
an[pk] = (long)(ap = apell(e, pl));
else
{
av = avma;
p1 = mulii(ap, an[pk / p]);
p2 = mulsi(p, an[pk / p / p]);
tetpil = avma;
an[pk] = lpile(av, tetpil, subii(p1, p2));
}
for(m = n / pk; m > 1; m--)
if (an[m] && (m % p)) an[m * pk] = lmulii(an[m], an[pk]);
}
}
return an;
}
GEN hell3(e,a,prec)
GEN e,a;
long prec;
{
long av=avma,tetpil;
GEN p1,p2,z,q,w,piisurw,pitemp;
if((typ(e)!=17)||(lg(e)<20)) err(elliper1);
if(!oncurve(e,a)) err(heller1);
if(lg(a)<3) return gzero;
constpi(prec);
z=zell(e,a,prec);pitemp=mppi(prec);piisurw=gdiv(gmul(gi,pitemp),w=(GEN)e[15]);
q=gexp(gmul(e[16],piisurw),prec);
p1=gdiv(theta(q,gmul(gdiv(z,w),pitemp),prec),thetanullk(q,1,prec));
p1=gmul(p1,gdiv(w,pitemp));
p1=glog(gmul(p1,gdiv(denom(a[1]),denom(a[2]))),prec);
p2=gimag(z);
if((gcmp0(p2))||((gexpo(p2)-gexpo(gimag(e[16])))<-2))
p1=gneg(p1);
else {p2=gmul(gmul2n(pitemp,-2),gimag(gdiv(e[16],e[15])));p1=gsub(p2,p1);}
tetpil=avma;return gerepile(av,tetpil,greal(p1));
}
GEN hell(e,a,prec)
GEN e,a;
long prec;
{
long av=avma,tetpil,n;
GEN p1,p2,y,z,q,psi2,pi2surw,pi2isurw,qn,ps;
if((typ(e)!=17)||(lg(e)<20)) err(elliper1);
pi2surw=gdiv(gmul2n(mppi(prec),1),(GEN)e[15]);
pi2isurw=cgetg(3,6);pi2isurw[1]=zero;pi2isurw[2]=(long)pi2surw;
z=gmul(greal(zell(e,a,prec)),pi2surw);
q=greal(gexp(gmul(e[16],pi2isurw),prec));
psi2=gadd(e[3],gadd(gmul(e[1],a[1]),gmul2n(a[2],1)));
y=gsin(z,prec);n=0;qn=gun;ps=gneg(q);
do
{
n++;p1=gsin(gmulsg(n+n+1,z),prec);qn=gmul(qn,ps);
ps=gmul(ps,q);p1=gmul(p1,qn);
y=gadd(y,p1);
}
while(gexpo(qn)>= -((prec-2)<<5));
p1=gmul(gsqr(gdiv(gmul2n(y,1),psi2)),pi2surw);
p2=gsqr(gsqr(gdiv(p1,gsqr(gsqr(denom(a[1]))))));
p1=gdiv(gmul(p2,q),e[12]);
p1=gmul2n(glog(gabs(p1,prec),prec),-5);
tetpil=avma;return gerepile(av,tetpil,gneg(p1));
}
GEN ghell(e,a,prec)
GEN e,a;
long prec;
{
long av=avma,av2,tetpil,lx,i,n,n2,grandn;
GEN p,p1,p2,x,y,z,phi2,psi2,psi3,logdep;
/* On suppose que la courbe est a coefficients entiers donne par un
modele minimal */
if((typ(e)!=17)||(lg(e)<20)) err(elliper1);
if(!oncurve(e,a)) err(heller1);
if(lg(a)<3) return gzero;
x=(GEN)a[1];y=(GEN)a[2];
psi2=numer(gadd(e[3],gadd(gmul(e[1],x),gmul2n(y,1))));
p2=gadd(gmulsg(3,e[7]),gmul(x,gadd(e[6],gmulsg(3,x))));
psi3=numer(gadd(e[9],gmul(x,gadd(gmulsg(3,e[8]),gmul(x,p2)))));
if((!signe(psi2))||(!signe(psi3))) {avma=av;return gzero;}
phi2=numer(gsub(gadd(e[4],gmul(x,gadd(shifti(e[2],1),gmulsg(3,x)))),gmul(e[1],
y)));
p1=(GEN)factor(mppgcd(psi2,phi2))[1];lx=lg(p1);
tetpil=avma;z=hell(e,a,prec);av2=avma;
if(lx==1) return gerepile(av,tetpil,z);
for(i=1;i<lx;i++)
{
p=(GEN)p1[i];n2=ggval(psi2,p);logdep=gneg(glog(p,prec));
if(signe(resii(e[10],p)))
{
grandn=ggval(e[12],p);n=n2<<1;
if(grandn)
{
if(n>grandn) n=grandn;
p2=divrs(mulsr(n*(grandn+grandn-n),logdep),(grandn<<3));
tetpil=avma;z=gadd(z,p2);
}
else avma=av2;
}
else
{
n=ggval(psi3,p);
if(n>=3*n2) p2=gdivgs(mulsr(n2,logdep),3);
else p2=gmul2n(mulsr(n,logdep),-3);
tetpil=avma;z=gadd(z,p2);
}
}
return gerepile(av,tetpil,z);
}
GEN ghell2(e,a,prec)
GEN e,a;
long prec;
{
long av=avma,tetpil,av2,lx,i,n,n2,grandn;
GEN p,p1,p2,x,y,z,phi2,psi2,psi3,logdep;
/* On suppose que la courbe est a coefficients entiers donne par un
modele minimal */
if((typ(e)!=17)||(lg(e)<20)) err(elliper1);
if(!oncurve(e,a)) err(heller1);
if(lg(a)<3) return gzero;
x=(GEN)a[1];y=(GEN)a[2];
psi2=numer(gadd(e[3],gadd(gmul(e[1],x),gmul2n(y,1))));
if(!signe(psi2)) return gzero;
phi2=numer(gsub(gadd(e[4],gmul(x,gadd(shifti(e[2],1),gmulsg(3,x)))),gmul(e[1],
y)));
p1=(GEN)factor(mppgcd(psi2,phi2))[1];lx=lg(p1);
tetpil=avma;z=hell2(e,a,prec);av2=avma;
if(lx==1) return gerepile(av,tetpil,z);
for(i=1;i<lx;i++)
{
p=(GEN)p1[i];n2=ggval(psi2,p);logdep=gneg(glog(p,prec));
if(signe(resii(e[10],p)))
{
grandn=ggval(e[12],p);n=n2<<1;
if(grandn)
{
if(n>grandn) n=grandn;
p2=divrs(mulsr(n*(grandn+grandn-n),logdep),(grandn<<3));
tetpil=avma;z=gadd(z,p2);
}
else avma=av2;
}
else
{
p2=gadd(gmulsg(3,e[7]),gmul(x,gadd(e[6],gmulsg(3,x))));
psi3=numer(gadd(e[9],gmul(x,gadd(gmulsg(3,e[8]),gmul(x,p2)))));
n=ggval(psi3,p);
if(n>=3*n2) p2=gdivgs(mulsr(n2,logdep),3);
else p2=gmul2n(mulsr(n,logdep),-3);
tetpil=avma;z=gadd(z,p2);
}
}
return gerepile(av,tetpil,z);
}
GEN ghell3(e,a,prec)
GEN e,a;
long prec;
{
long av=avma,tetpil,av2,lx,i,n,n2,grandn;
GEN p,p1,p2,x,y,z,phi2,psi2,psi3,logdep;
/* On suppose que la courbe est a coefficients entiers donne par un
modele minimal */
if((typ(e)!=17)||(lg(e)<20)) err(elliper1);
if(!oncurve(e,a)) err(heller1);
if(lg(a)<3) return gzero;
x=(GEN)a[1];y=(GEN)a[2];
psi2=numer(gadd(e[3],gadd(gmul(e[1],x),gmul2n(y,1))));
if(!signe(psi2)) return gzero;
phi2=numer(gsub(gadd(e[4],gmul(x,gadd(shifti(e[2],1),gmulsg(3,x)))),gmul(e[1],
y)));
p1=(GEN)factor(mppgcd(psi2,phi2))[1];lx=lg(p1);
tetpil=avma;z=hell3(e,a,prec);av2=avma;
if(lx==1) return gerepile(av,tetpil,z);
for(i=1;i<lx;i++)
{
p=(GEN)p1[i];n2=ggval(psi2,p);logdep=gneg(glog(p,prec));
if(signe(resii(e[10],p)))
{
grandn=ggval(e[12],p);n=n2<<1;
if(grandn)
{
if(n>grandn) n=grandn;
p2=divrs(mulsr(n*(grandn+grandn-n),logdep),(grandn<<3));
tetpil=avma;z=gadd(z,p2);
}
else avma=av2;
}
else
{
p2=gadd(gmulsg(3,e[7]),gmul(x,gadd(e[6],gmulsg(3,x))));
psi3=numer(gadd(e[9],gmul(x,gadd(gmulsg(3,e[8]),gmul(x,p2)))));
n=ggval(psi3,p);
if(n>=3*n2) p2=gdivgs(mulsr(n2,logdep),3);
else p2=gmul2n(mulsr(n,logdep),-3);
tetpil=avma;z=gadd(z,p2);
}
}
return gerepile(av,tetpil,z);
}
GEN lseriesell(e,s,N,A,prec)
GEN e,s,N,A;
long prec;
{
long av=avma,tetpil,l,n,eps;
GEN z,p1,p2,cg,cg1,v,cga,cgb,s2,ns,gs;
if((typ(e)!=17)||(lg(e)<14)||(typ(N)!=1)||(!signe(N))) err(elliper1);
if(gsigne(A)<=0) err(lserieser1);
if(gcmpgs(A,1)<0) A=ginv(A);
eps=signe(N);if(eps<0) N=gneg(N);
cg1=mppi(prec);setexpo(cg1,2);cg=divrr(cg1,gsqrt(N,prec));
cga=gmul(cg,A);cgb=gdiv(cg,A);
l=(C2*(prec-2)+fabs(gtodouble(s)-1.)*log(rtodbl(cga)))/rtodbl(cgb)+1;
v=anell(e,l);
s2=gsubsg(2,s);ns=gpui(cg,gsubgs(gmul2n(s,1),2),prec);
z=gzero;
if(typ(s)==1)
{
if(signe(s)>0) gs=mpfactr(itos(s)-1,prec);
else {avma=av;return gzero;}
}
else gs=ggamma(s,prec);/* gs2=ggamma(s2,prec); */
for(n=1;n<=l;n++)
{
p1=gdiv(incgam4(s,gmulsg(n,cga),gs,prec),gpui(stoi(n),s,prec));
p2=gdiv(gmul(ns,incgam(s2,gmulsg(n,cgb),prec)),gpui(stoi(n),s2,prec));
if(eps<0) p2=gneg(p2);
z=gadd(z,gmul(gadd(p1,p2),v[n]));
}
tetpil=avma;return gerepile(av,tetpil,gdiv(z,gs));
}
/********************************************************************/
/********************************************************************/
/** **/
/** Algorithme de Tate **/
/** (cf Anvers IV) **/
/** Type de Kodaira, modele minimal global **/
/** **/
/********************************************************************/
/********************************************************************/
/*
Etant donnes une courbe elliptique sous forme longue e, dont les coefficients
sont entiers, et un nombre premier p1, renvoie le type de la fibre en p du
modele de Neron de la courbe, ainsi que les changements de variables
necessaires, sous la forme d'un triplet [f, kod, v].
L'entier f est l'exposant du conducteur.
l'entier kod, est le type de kodaira. les types II, III et IV sont codes 2,
3, et 4 respectivement. Les types II*, III* et IV* donnent -2, -3 et -4. Le
type I0* donne -1, les types Inu et Inu* donnent 4+nu et -4-nu. Enfin, le
type I0 donne 1.
v est un quadruplet [u, r, s, t] qui permet de passer a un modele minimal.
En general, on ne s'interessera a ce vecteur que si u <> 1.
L'algorithme est bien sur celui de Tate dans Anvers IV. Compte tenu des
remarques du bas de la page 46, l'algorithme long n'est utilise que pour
p = 2 ou p = 3.
Remarque. L'algorithme donne aussi le nombre c (voir Tate), qui n'est pas
calcule ici. Peut-etre faudrait-il le faire ?
*/
void cumule(), cumule1();
GEN localreduction(e, p1)
GEN e, p1;
{
long av = avma, tetpil, k, p, f, kod, nu, nuj, nudelta;
int a21, a42, a63, a32, a64, theroot;
GEN pk, p2k, pk1, p4, p6, a2prime, a3prime;
GEN p2, p3, r = gzero, s = gzero, t = gzero, v, result;
if ((typ(e) != 17) || (lg(e) < 14)) err(elliper1);
v = cgetg(5,17); v[1] = un; v[2] = v[3] = v[4] = zero;
nudelta = ggval(e[12], p1);
if (gcmpgs(p1, 3) > 0) /* p different de 2 ou 3 */
{
nuj = gcmp0(e[13]) ? 0 : - ggval(e[13], p1);
k = (nuj > 0 ? nudelta - nuj : nudelta) / 12;
if (k > 0) /* modele non minimal */
{
pk = gpuigs(p1, k);
if (mpodd(e[1]))
s = shifti(subii(pk, e[1]), -1);
else
s = negi(shifti(e[1], -1));
p2k = mulii(pk, pk);
a2prime = subii(e[2], mulii(s, addii(e[1], s)));
switch(itos(modis(a2prime, 3)))
{
case 0: r = negi(divis(a2prime, 3)); break;
case 1: r = divis(subii(p2k, a2prime), 3); break;
case 2: r = negi(divis(addii(a2prime, p2k), 3)); break;
default: err(tater1);
}
a3prime = addii(e[3], mulii(r, e[1]));
if (mpodd(a3prime))
t = shifti(subii(mulii(pk, p2k), a3prime), -1);
else
t = negi(shifti(a3prime, -1));
v[1] = (long)pk; v[2] = (long)r; v[3] = (long)s; v[4] = (long)t;
nudelta -= 12 * k;
}
if (nuj > 0) switch(nudelta - nuj)
{
case 0: f = 1; kod = 4 + nuj; break; /* Inu */
case 6: f = 2; kod = - 4 - nuj; break; /* Inu* */
default: err(tater1);
}
else switch(nudelta)
{
case 0: f = 0; kod = 1; break; /* I0, regulier */
case 2: f = 2; kod = 2; break; /* II */
case 3: f = 2; kod = 3; break; /* III */
case 4: f = 2; kod = 4; break; /* IV */
case 6: f = 2; kod = -1; break; /* I0* */
case 8: f = 2; kod = -4; break; /* IV* */
case 9: f = 2; kod = -3; break; /* III* */
case 10: f = 2; kod = -2; break; /* II* */
default: err(tater1);
}
}
else for(;;)
{
/* ici, p = 2 ou p = 3 */
if (!nudelta) {f = 0; kod = 1; goto exit;}
/* I0 */
if (!divise(e[6], p1)) {f = 1; kod = 4 + nudelta; goto exit;}
/* Inu */
p = itos(p1);
if (p == 2)
{
r = modis(e[4], 2);
s = modis(addii(r, e[2]), 2);
if (signe(r)) t = modis(addii(addii(e[4], e[5]), s), 2);
else t = modis(e[5], 2);
}
else /* p == 3 */
{
r = negi(modis(e[8], 3));
s = modis(e[1], 3);
t = modis(addii(e[3], mulii(e[1], r)), 3);
}
cumule(&v, &e, gun, r, s, t); /* p | a1, a2, a3, a4 et a6 */
p2 = stoi(p*p);
if(!divise(e[5], p2)) {f = nudelta; kod = 2; goto exit;}
/* II */
p3 = stoi(p*p*p);
if(!divise(e[9], p3)) {f = nudelta - 1; kod = 3; goto exit;}
/* III */
if (!divise(e[8], p3)) {f = nudelta - 2; kod = 4; goto exit;}
/* IV */
/* now for the last five cases... */
if (!divise(e[5], p3))
cumule(&v, &e, gun, gzero, gzero, p == 2 ? stoi(2) : modis(e[3], 9));
/* p | a1, a2; p^2 | a3, a4; p^3 | a6 */
a21 = aux((GEN)e[2], p, 1); a42 = aux((GEN)e[4], p, 2);
a63 = aux((GEN)e[5], p, 3);
switch (numroots3(a21, a42, a63, p, &theroot))
{
case 3: f = nudelta - 4; kod = -1; goto exit;
/* I0* */
case 2: /* calcul de nu */
if (theroot) cumule(&v, &e, gun, stoi(theroot * p), gzero, gzero);
/* p | a1; p^2 | a2, a3; p^3 | a4; p^4 | a6 */
nu = 1;
pk = p2;
p2k = stoi(p * p * p * p);
for(;;)
{
if (numroots2(1, aux2((GEN)e[3], p, pk),
-aux2((GEN)e[5], p, p2k), p, &theroot) == 2)
break;
if (theroot) cumule(&v, &e, gun, gzero, gzero, mulsi(theroot,pk));
pk1 = pk; pk = mulsi(p, pk); p2k = mulsi(p, p2k);
nu++;
if (numroots2(a21, aux2((GEN)e[4], p, pk),
aux2((GEN)e[5], p, p2k), p, &theroot) == 2)
break;
if (theroot) cumule(&v, &e, gun, mulsi(theroot, pk1), gzero, gzero);
p2k = mulsi(p, p2k);
nu++;
}
f = nudelta - 4 - nu; kod = -4 - nu; goto exit;
/* Inu* */
case 1:
if (theroot) cumule(&v, &e, gun, stoi(theroot * p), gzero, gzero);
/* p | a1; p^2 | a2, a3; p^3 | a4; p^4 | a6 */
a32 = aux((GEN)e[3], p, 2); a64 = aux((GEN)e[5], p, 4);
if(numroots2(1, a32, -a64, p, &theroot) == 2)
{f = nudelta - 6; kod = -4; goto exit;}
/* IV* */
if (theroot) cumule(&v, &e, gun, gzero, gzero, stoi(theroot*p*p));
/* p | a1; p^2 | a2; p^3 | a3, a4; p^5 | a6 */
p4 = mulii(p2, p2);
if (!divise(e[4], p4)) {f = nudelta - 7; kod = -3; goto exit;}
/* III* */
p6 = mulii(p4, p2);
if (!divise(e[5], p6)) {f = nudelta - 8; kod = -2; goto exit;}
/* II* */
cumule(&v, &e, p1, gzero, gzero, gzero); /* non minimal, on repart
pour un tour */
nudelta -= 12;
}
}
exit:
tetpil = avma;
result = cgetg(4, 17);
result[1] = lstoi(f); result[2] = lstoi(kod); result[3] = lcopy(v);
return gerepile(av, tetpil, result);
}
/*
Calcul de toutes les fibres non elliptiques d'une courbe sur Z.
Etant donne une courbe elliptique sous forme longue e, dont les coefficients
sont entiers, renvoie une matrice dont les lignes sont de la forme
[p, fp, kodp]. Il y a une ligne par diviseur premier du discriminant.
Ceci n'est pas implemente dans Gp. Est-ce utile ?
*/
GEN globaltatealgo(e)
GEN e;
{
long k, l,av;
GEN p1, p2, p3, prims, result;
if ((typ(e) != 17) || (lg(e) < 14)) err(elliper1);
prims = decomp(e[12]);
l = lg(p1 = (GEN)prims[1]);
p2 = (GEN)prims[2];
if ((long)prims == avma) cgiv(prims);
result = cgetg(4, 19);
result[1] = (long)p1;
result[2] = (long)p2;
result[3] = (long)(p3 = cgetg(l, 18));
for(k = 1; k < l; k++) p3[k] = lgeti(3);
av = avma;
for(k = 1; k < l; k++)
{
GEN q = localreduction(e, (GEN)p1[k]);
affii(q[1], p2[k]);
affii(q[2], p3[k]);
avma = av;
}
return result;
}
/*
Algorithme de reduction d'une courbe sur Q a sa forme standard.
Etant donne une courbe elliptique sous forme longue e, dont les coefficients
sont rationnels, renvoie son [N, [u, r, s, t]], ou N est le conducteur
arithmetique de e, et [u, r, s, t] est le changement de variables qui reduit
e a sa forme minimale globale dans laquelle a1 et a3 valent 0 ou 1, et a2
vaut -1, 0 ou 1 et tel que u est un rationnel positif.
*/
GEN globalreduction(e1)
GEN e1;
{
long i, k, l, m, tetpil, av = avma;
GEN p1, prims, result, N = gun, u = gun, r, s, t;
GEN v = cgetg(5, 17), a = cgetg(7, 17), e = cgetg(20, 17);
if ((typ(e1) != 17) || (lg(e1) < 14)) err(elliper1);
for(i = 1; i < 5; i++) a[i] = e1[i]; a[5] = zero; a[6] = e1[5];
prims = decomp(denom(a));
l = lg(p1 = (GEN)prims[1]);
for(k = 1; k < l; k++)
{
int n = 0;
for(i = 1; i < 7; i++)
if (!gcmp0(a[i]))
{
m = i * n + ggval(a[i], (GEN)p1[k]);
while(m < 0) {n++; m += i;}
}
u = gmul(u, gpuigs(p1[k], n));
}
v[1] = linv(u); v[2] = v[3] = v[4] = zero;
for(i = 1; i < 14; i++) e[i] = e1[i];
for(; i < 20; i++) e[i] = zero;
if (!gcmp1(u)) e = coordch(e, v);
prims = decomp(e[12]);
l = lg(p1 = (GEN)prims[1]);
for(k = (signe(e[12]) < 0) + 1; k < l; k++)
{
GEN q = localreduction(e, (GEN)p1[k]);
GEN v1 = (GEN)q[3];
N = mulii(N, gpui(p1[k], q[1]));
if (!gcmp1((GEN)v1[1])) cumule1(&v, &e, v1);
}
s = gdiventgs(e[1], -2);
r = gdiventgs(gaddgs(gsub(gsub(e[2], gmul(s, e[1])), gsqr(s)), 1), -3);
t = gdiventgs(gadd(e[3], gmul(r, e[1])), -2);
cumule(&v, &e, gun, r, s, t);
tetpil = avma;
result = cgetg(3, 17); result[1] = lcopy(N); result[2] = lcopy(v);
return gerepile(av, tetpil, result);
}
/* renvoie a_{k,l} avec les notations de Tate */
static int aux(ak, p, l)
GEN ak;
int p, l;
{
long av = avma, pl = p, res;
while(--l) pl *= p;
res = itos(modis(divis(ak, pl), p));
avma = av;
return res;
}
static int aux2(ak, p, pl)
GEN ak, pl;
int p;
{
long av = avma, res;
res = itos(modis(divii(ak, pl), p));
avma = av;
return res;
}
/* renvoie le nombre de racines distinctes du polynome XXX + aXX + bX + c
modulo p s'il y a une racine multiple, elle est renvoyee dans *mult */
static int numroots3(a, b, c, p, mult)
int a, b, c, p, *mult;
{
if (p == 2)
if ((c + a * b) & 1) return 3;
else {*mult = b; return (a + b) & 1 ? 2 : 1;}
else
if (a % 3) {*mult = a * b; return (a * b * (1 - b) + c) % 3 ? 3 : 2;}
else {*mult = -c; return b % 3 ? 3 : 1;}
}
/* idem pour aXX +bX + c */
static int numroots2(a, b, c, p, mult)
int a, b, c, p, *mult;
{
if (p == 2) {*mult = c; return b & 1 ? 2 : 1;}
else {*mult = a * b; return (b * b - a * c) % 3 ? 2 : 1;}
}
/* cumule les effets de plusieurs chgts de variable. On traite a part les cas
particuliers frequents, tels que soit u = 1, soit r' = s' = t' = 0 */
static void cumulev(vtotal, u, r, s, t)
GEN *vtotal, u, r, s, t;
{
long av = avma, tetpil;
GEN temp, v = *vtotal, v3 = cgetg(5, 17);
if (gcmp1(v[1]))
{
v3[1] = lcopy(u);
v3[2] = ladd(v[2], r);
v3[3] = ladd(v[3], s);
av = avma;
temp = gadd(v[4], gmul(v[3], r));
tetpil = avma;
v3[4] = lpile(av, tetpil, gadd(temp, t));
}
else if (gcmp0(r) && gcmp0(s) && gcmp0(t))
{
v3[1] = lmul(v[1], u);
v3[2] = lcopy(v[2]);
v3[3] = lcopy(v[3]);
v3[4] = lcopy(v[4]);
}
else /* cas general */
{
v3[1] = lmul(v[1], u);
temp = gmul(v[1],v[1]);
v3[2] = ladd(gmul(temp, r), v[2]);
v3[3] = ladd(gmul(v[1], s), v[3]);
v3[4] = ladd(v[4], gmul(temp, gadd(gmul(v[1], t), gmul(v[3], r))));
tetpil = avma;
v3 = gerepile(av, tetpil, gcopy(v3));
}
*vtotal = v3;
}
static void cumule(vtotal, e, u, r, s, t)
GEN *vtotal, *e, u, r, s, t;
{
long av = avma, tetpil;
GEN v2 = cgetg(5, 17);
v2[1] = (long)u; v2[2] = (long)r; v2[3] = (long)s; v2[4] = (long)t;
tetpil = avma;
*e = gerepile(av, tetpil, coordch(*e, v2));
cumulev(vtotal, u, r, s, t);
}
static void cumule1(vtotal, e, v2)
GEN *vtotal, *e, v2;
{
*e = coordch(*e, v2);
cumulev(vtotal, (GEN)v2[1], (GEN)v2[2], (GEN)v2[3], (GEN)v2[4]);
}